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ABSTRACT 

We present a determination of the bispectrum of the flux in the Lyman-cv for- 
est of QSO absorption spectra obtained from LUQAS which consists of spectra 
observed with the high resolution Echelle spectrograph UVES. Typical errors 
on the observed bispectrum as obtained from a jack-knife estimator are ~ 50 
%. For wavenumbers in the range 0.03 s/km < k < 0.1 s/km the observed bis- 
pectrum agrees within the errors with that of the synthetic absorption spectra 
obtained from numerical hydro-simulations of a ACDM model with and with- 
out feedback from star formation. Including galactic feedback changes the 
bispectrum by less than 10%. At smaller wavenumbers the associated metal 
absorption lines contribute about 50% to the bispectrum and the observed 
bispectrum exceeds that of the simulations. At wavenumbers k < 0.03 s/km 
second-order perturbation theory applied to the flux spectrum gives a rea- 
sonable (errors smaller than 30%) approximation to the bispectra of observed 
and simulated absorption spectra. The bispectrum of the observed absorp- 
tion spectra also agrees, within the errors, with that of a randomized set of 
absorption spectra where a random shift in wavelength has been added to ab- 
sorption lines identified with VPFIT. This suggests that for a sample of the 
size presented here, the errors on the bispectrum are too large to discriminate 
between models with very different 3D distribution of Lyman-a absorption. If 
it were possible to substantially reduce these errors for larger samples of ab- 
sorption spectra, the bispectrum might become an important statistical tool 
for probing the growth of gravitational structure in the Universe at redshift 
z> 2. 

Key words: Cosmology: intergalactic medium - large-scale structure of uni- 
verse - quasars: absorption lines 



1 INTRODUCTION 

The Lyman-a forest in the absorption spectra of high- 
redshift QSOs has been recognized as an important 
tool for studying the underlying matter distribution (see 

* Based on data taken from the ESO archive obtained with 
UVES at VLT, Paranal, Chile. 



Rauch 1998 and Weinberg et al. 1999 for excellent re- 
views). In a previous paper (Kim et al. 2003, K03) an 
analysis of the flux power-spectrum of a large sample of 
UVES QSO absorption spectra (LUQAS) has been pre- 
sented. Here we investigate the bispectrum of the flux 
in a subset of these absorption spectra. 

The bispectrum is the Fourier transform of the 
three-point correlation function. As the flux in the 
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Lyman-a forest is a sensitive probe of the matter dis- 
tribution it should probe the topology of the 3D matter 
distribution in more detail than the two-point corre- 
lation function or the power-spectrum. It also can be 
used as a complementary tool to determine cosmologi- 
cal parameters (e.g. Fry 1994; Verde et al. 2002) and 
maybe also the physical state of the IGM. Gravita- 
tional growth induces correlations between large scale 
modes and small scale power which can be probed by the 
bispectrum. Zaldarriaga et al. (2001) pointed out that 
these correlations may be used to discriminate between 
fluctuations due to large-scale structure in the matter 
distribution and those produced by non-gravitational 
processes such as fluctuations in the continuum emission 
of the quasar. Mandelbaum et al. (2003) showed that 
with the SDSS (Sloan Digital Sky Survey) QSO sample 
it should be possible to use higher order statistics such 
as the bispectrum to determine amplitude, slope and 
curvature of the slope of the matter power spectrum 
with an accuracy of a few percent provided systematic 
errors are under control. 

The plan of the paper is as follows. In Section 2 
we briefly describe our sample. In Section 3 we de- 
fine the flux bispectrum and compare the data with 
hydro-dynamical simulations, an analytical prediction 
and randomized spectra. Section 4 contains a summary 
of our results. 
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Figure 1. Filled circles show the one-dimensional bispec- 
trum of the Lyman-a forest flux of the LUQAS sample in the 
rcdshift range 2 < z < 2.4. The corresponding bispectrum 
of the associated metal absorption is shown by filled dia- 
monds. Also shown is the flux bispectrum of synthetic spec- 
tra computed from hydro-dynamical simulations of a ACDM 
model at z = 2.3 (see Section 3.3). Filled triangles are for 
a simulation with a box size of 33.75 Mpc which in- 

clude feedback from galactic winds. Empty triangles are for 
a simulation with box size of 10 h~ 1 Mpc without galactic 
feedback. Squares indicate positive values of B(k). Errors are 
obtained with a jack-knife estimate. 



2 THE DATA 

The LUQAS sample consists of 27 spectra taken with 
the Ultra- Violet Echelle Spectrograph (UVES) on VLT 
(Paranal, Chile) over the period 1999-2002. The spec- 
tra were drawn from the ESO archive and are publicly 
available to the ESO community. The total redshift path 
of the sample is Az = 13.75. The median redshift of the 
sample is (z) = 2.25 and the number of spectra covering 
the median redshift is 17. The analysis of the bispectrum 
has been restrained to the redshift range 2 < z < 2.4, 
where approximately 15 QSOs contribute. The reason 
for this is the strong dependence of the bispectrum on 
the mean flux. The typical S/N is ~ 50. 
The 

spectra were reduced with the ECHELLE/UVES en- 
vironment of the software package MIDAS. Lists of 
HI absorption lines and associated metal absorption 
lines identified and fitted with VPFIT (Carswell et ah: 
http://www.ast.cam.ac.uk/~rfc/vpfit.html) were avail- 
able for 8 and 13 QSOs, respectively. Five out of the 
eight QSOs with lists of HI absorption are in the red- 
shift range for which we calculated the bispectrum of 
the flux. We have removed four damped/sub-damped 
Lyman-a systems before computing the bispectrum. For 
a more complete description of the sample and for the 
analysis of the flux power-spectrum and its evolution 
with redshift we refer to K03. 



3 THE BISPECTRUM OF THE 
Lyman-a FOREST 

3.1 Method 

We consider the quantity 8f — (F — (F))/{F), where 
F is the continuum fitted transmitted flux and (F) the 
average flux for each QSO (Hui et al. 2001, Croft et al. 
2002). Note that 5f corresponds to the flux estimator 
F 2 of K03. 

We follow the definition of Matarrese, Verde & 
Heavens (1997) and Verde et al. (1998) for the bispec- 
trum. We use the real part of the three point function 
in A;— space , Df = Re(<5F(fci) (^(fe) (Jf^)), for closed 
triangles ki + ki + ks = 0. Sp(k) is the Fourier trans- 
form of SF. Df is related to the bispectrum of the flux 
B F (k 1 ,k 2 , k s ) 

(D F ) = 2itB F (k 1 ,k2,k 3 )S D (k 1 + k2 + k 3 ) . (1) 

5 D (k) is the one-dimensional Dirac delta function and 
(■) indicates the ensemble average. We compute the one- 
dimensional bispectrum. Our triangles are thus degen- 
erate and we choose the configuration for which ki = fc2 
and &3 = —2k\. In the following we will always show 
the flux bispectrum as a function of the wavenumber 
k — k±. 
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3.2 The observed bispectrum 

The filled cirles in Fig. 1 show the one-dimensional bis- 
pectrum of the Lyman-a forest flux of the LUQAS sam- 
ple in the redshift range 2 < z < 2.4 plotted in di- 
mensionless units (Table 1). The bin size was chosen 
such that visual comparison of different bispectra in the 
plots is not hindered by too large bin-to-bin fluctua- 
tions. As in K03 the error bars are computed with a 
jack-knife estimator, which gives results very similar to 
a bootstrap resampling of the data. Note that the bis- 
pectrum computed from the continuum fitted observed 
absorption spectra is negative. The same was found by 
Mandelbaum et al. (2003) for a statistic similar to the 
bispectrum in their analysis of synthetic spectra from 
numerical simulations. The negative sign of the bispec- 
trum is expected if the higher order correlations arise 
from gravitational growth (Zaldarriaga et al. 2001). 

For 13 of the 27 spectra we have metal line lists. 
From these lists we have produced artificial spectra 
which contain the metal lines only. The diamonds in 
Figure 1 show the mean bispectrum of these metal-line- 
only spectra which have a median redshift z — 2.36. For 
k > 0.1 s/km the identified metal lines contribute up 
to 50% to the bispectrum. At larger scales their contri- 
bution is considerably smaller. 

We have also computed the flux bispectrum for the 
flux estimator F3 of K03 (see Figure 2 of K03) , which is 
defined as 8f3 = F ' /(F) — 1, where F is the flux of the 
spectrum without continuum fitting. At k < k con t ~ 
0.003 s/km continuum fluctuations dominate the flux 
bispectrum and the values of the bispectrum become 
positive. The same was found by K03 for the flux power- 
spectrum. We therefore do not plot the bispectrum for 
k < 0.003 s/km. 

3.3 Comparison with hydro-dynamical 
simulations 

We have also calculated the bispectrum for synthetic 
spectra obtained from two outputs at z = 2.3 from 
the large set of hydro-dynamical simulations presented 
by Springel & Hernquist (2003a) . The simulations were 
performed with the code GADGET (Springel, Yoshida 
& White 2001) modified to conserve entropy (Springel 
& Hernquist 2002) and are for a ACDM model with 
parameters () = 0.3 and $1a = 0.7, Hubble constant 
Ho = 100/ikms -1 Mpc -1 with h = 0.7, baryon density 
fib = 0.04, a power spectrum with primordial spectral 
index n = 1, and rras fluctuation amplitude on 
Mpc scale, a s = 0.9. 

The simulation outputs are part of the 03 run and 
the D4 run. The 03 simulation has a box size of 10 
Mpc box size, no feedback due to star formation and 
was performed with 2 x 144 particles. The D4 run has 
a box size of 33.75 h^ 1 Mpc, 2 x 216 3 particles and in- 
cludes feedback by galactic winds. The mass of a gas 
particle is 3.72 x 10 6 M Q and 4.24 x 1O 7 M in 03 and 
D4, respectively. Both simulations assume an UV back- 
ground as modelled by Haardt & Madau (1996) and 
follow star formation with a hybrid multi-phase model 
(see Springel & Hernquist 2003b for details). Synthetic 



spectra are extracted from the simulation in the usual 
way and are normalized to reproduce the effective op- 
tical depth T eff = -ln(F) = 0.165 at (z) = 2.2 in the 
observed spectra used for calculating the bispectrum. 

The triangles in Fig. 1 show the flux bispectrum of 
the synthetic spectra with and without feedback from 
galactic winds. At k < 0.1 s/km the difference between 
the two simulations is smaller than 10 %. This confirms 
that the feedback from high-redshift galaxies is expected 
to have little effect on the gas responsible for Lyman- 
a absorption (Theuns et al. 2002). This is most likely 
because the simulations mainly differ in the distribution 
of hot gas which has a small filling factor even for the 
simulation with feedback. The difference in box size and 
resolution also appear to have a little effect. For k < 0.1 
s/km the bispectrum of synthetic and observed spectra 
agree within the errors. At larger wavenumbers the ob- 
served bispectrum exceeds that of the synthetic spectra. 
This is expected since at these wavenumbers the metal 
lines contribute significantly to the bispectrum (see Fig- 
ure 1). 



3.4 The bispectrum for randomized 
absorption spectra 

Viel et al. (2003b) have demonstrated that the Lyman- 
a forest flux power spectrum of "randomised" QSO ab- 
sorption spectra is comparable in shape and amplitude 
to the flux power spectrum of the original observed spec- 
tra. They found that this is because of a large contri- 
bution of the shape of the Voigt profiles of discrete ab- 
sorption systems to the flux power spectrum. We have 
performed here a similar test for the bispectrum. For 5 
of the observed spectra in the redshift range where we 
have calculated the flux bispectrum we have a list of HI 
absorption lines obtained from Voigt profile fitting with 
VPFIT. For these we have produced 50 artificial ran- 
domized spectra each, where we have randomly shifted 
the hydrogen lines in wavelength. In Fig. 2 we compare 
the corresponding flux bispectrum to the flux bispec- 
trum of the observed spectra. At k < 0.1 s/km they 
agree within the errors. A random superposition of Voigt 
profiles appears to reproduce the observed bispectrum 
similarly well as the numerical simulations. Metal lines 
were omitted in the reshuffled spectra hence the dis- 
crepancy at large wavenumbers. As in Viel et al. (2003) 
we have also calculated the bispectrum for randomized 
spectra where we halved and doubled the Doppler pa- 
rameter of the lines. The result is similar. Broadening 
and narrowing of the lines leads to an approximately lin- 
ear shift in wavenumber. This suggests that the bispec- 
trum measures to a significant extent the shape of Voigt 
profiles. Note that the errors of the flux bispectrum are 
larger than that of the flux power spectrum. It is thus 
more difficult to discriminate between the flux bispec- 
trum of randomized spectra and the flux distribution of 
a CDM-like model than between the corresponding flux 
power spectra. 
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Figure 2. Comparison between the flux bispectrum of ob- 
served (filled circles) and "randomized" absorption spectra 
(empty triangles). The randomized absorption spectra have 
been obtained by randomly shifting the position of the ab- 
sorption lines as obtained with VPFIT within the observed 
wavelength range. Error bars are jack-knife estimates. 



3.5 A second order perturbation theory 
approximation 

The structures which give rise to the Lyman-a forest 
have been shown to be mildly non-linear and to be rea- 
sonably well described by simple approximation schemes 
such as the log-normal model, at least on large scales 
(Bi & Davidsen 1997; Viel et al. 2002a, 2002b, 2003a; 
Matarrese & Mohayaee 2002). This motivated the use 
of perturbation theory to obtain analytical approxima- 
tions in the Fluctuating Gunn-Peterson Approximation 
to get some insight into the shape and scaling of the the 
flux power spectrum at large scales (e.g. Zaldarriaga, 
Scoccimarro & Hui 2002). In the same spirit we will 
here make an attempt to use second-order perturbation 
theory (e.g. Matarrese et al. 1997) to derive an analyt- 
ical expression which relates flux power spectrum and 
flux bispectrum and test it with synthetic spectra. 

If the initial fluctuations are Gaussian and struc- 
tures grow via gravitational instability the three-point 
correlation function is a second-order quantity. I In the 
Fluctuating Gunn-Peterson Approximation (FGPA; e.g. 
Hernquist et al. 1996) the flux can be related to the den- 
sity field as 



exp[-A{l + 5 IG Mf] 



(2) 



where A and f3 = 2 — 0.7(7 — 1) depend on z. 7 is 
the power-law index of the gas temperature-density re- 
lation as usual. This assumes that redshift space dis- 
tortions are not important (but see the extra effect of 
this assumption in Verde et al. 1998) and neglects ther- 
mal broadening and instrumental noise. Moreover, we 
assume that at the scales of interest Sigm ~ &dm ~ 5. 



A more refined treatment, where these approximations 
are dropped, will be presented elsewhere. 

Expanding the density field to second order as 
<5(x) » <5 (1) (x) + <5 (2) (x), and using the FGPA we get 



^«6it«W(x)+«W(x)] + ^ a (x), 



(3) 



with 61 = -Aj3andb 2 = —A (3 (/3—1—A 0). The second- 
order quantity <5' 2 ' can be expressed as a quadratic 
combination of linear perturbations <5 < - 1 ' with a suit- 
able convolution kernel. The kernel is obtained by an 
expansion of the equations of gravitational instability 
up to second order (e.g. Catelan et al. 1995 and refer- 
ences therein). The expression for the Lyman-a bispec- 
trum follows then by projecting the 3D bispectrum (e.g. 
Matarrese et al. 1997) along the line-of-sight, 

B(ki,k 2 ,k 3 ) = ci + p{k 1 )p(k 2 ) 

+ c 1 [(k 1 k 2 -h 2 1 )p i - 1) (k 1 )p(k 2 ) 

+ {k2kx-"^kl)p ( - 1 \k2)p{k 1 ) 

+ ktp^(k 1 ) P ^ 1 Hk2)} + cyc.(l, 2, 

where ci = l/&i, c 2 — b 2 jb\. p(k) is the ID flux power 
spectrum. The spectral moment p^ _1 '(fc) is given by, 

/•oo 

p {e) {k) = \k\ 2l p{k)+2l dqq-^piq) (5) 

J\k\ 

with t = — 1. Note that in the above eqs. ki,k 2 , k$ have 
to be taken including their signs. The theoretical rela- 
tion between the flux bispectrum and power-spectrum 
only depends upon the IGM parameters A and j3. The 
very weak dependence on the matter density parameter 
Q m (e.g. Bouchet et al. 1992; Catelan et al. 1995) has 
been neglected as usual. 

In Figure 3 we compare the bispectrum calculated 
from the ID flux power spectrum of the sythentic spec- 
tra using eqs. (4) and (5) with the flux bispectrum cal- 
culated directly from the synthetic spectra. The rea- 
son for choosing the simulated and not the observed 
spectra to test the accuracy of the analytical relation 
between bispectrum and power-spectrum is the signifi- 
cantly smaller errors of the bispectrum of the synthetic 
spectra. We have used synthetic spectra calculated from 
the D4 run of the numerical simulations(see Section 
3.3). At k < 0.03 s/km the relation between bispectrum 
and power spectrum based on second-order perturba- 
tion theory is in good agreement with the simulations if 
we choose A = 0.48. The dashed and continuous curve 
are for 7 = 1.1 and 7 = 1.5, respectively. The differ- 
ences between a value of 7 = 1.1 and 7 = 1.5 are very 
small. At k > 0.03 s/km the analytical prediction and 
the true bispectrum start to differ. This is probably due 
to neglect of the effects of thermal broadening and Jeans 
smoothing. 

Note that in order to reproduce the correct ampli- 
tude we had to choose a different value A in the ana- 
lytical approximation of the bispectrum than that used 
in the synthetic spectra where A was set to reproduce 
the observed mean flux decrement. The synthetic spec- 
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Figure 3. The flux bispectrum of the Lyman-a forest calcu- 
lated from the flux power spectrum of synthetic spectra using 
the analytical relation of eq. (4). The continuous, dashed and 
dotted curves are for different parameter A and 7 of the Fluc- 
tuating Gunn- Peterson Approximation (see text for details). 
Triangles show the bispectrum calculated directly from the 
synthetic spectra. 

tra were calculated with A — 0.36 and the analytical 
prediction of the bispectrum for this value is shown by 
the dotted curve in Figure 3. There is a clear discrep- 
ancy in the amplitude. We have tried to find the reason 
why the analytical prediction gives a good match for a 
somewhat different value of A without much success. We 
can only speculate that non-linear effects, not properly 
taken into account by second-order perturbation theory, 
are responsible for that. 



4 RESULTS AND DISCUSSION 

We have computed the flux bispectrum of the LUQAS, 
sample a set of high resolution quasar spectra taken with 
UVES in the redshift range 2 < z < 2.4, and compared 
to the flux bispectrum of synthetic spectra of a hydro- 
dynamical simulation and that of randomized spectra. 
The main results are as follows: 

• at wavenumbers k < 0.1 s/km the flux bispectrum 
of the observed spectra and that of synthetic spectra 
obtained of hydro-dynamical simulations of a ACDM 
model agree well within the errors (50%); 

• including feedback from galactic winds has little 
effect on the flux bispectrum; 

• at wavenumbers k > 0.1 s/km identified metal lines 
contribute significantly to the flux bispectrum (of the 
order 50 %), while at larger scales the contribution by 
metal lines is negligible; 

• at wavenumbers k < 0.1 s/km the observed flux 
bispectrum is consistent within the errors with that ob- 
tained for randomized spectra; 



Table 1. Flux bispectrum( a ) of the LUQAS sample in the 
redshift range 2.0 < z < 2.4. 



k (s/km) 


B(k) (km/s) 2 


0.0030 


-2845.669 ± 1057.042 


0.0057 


-1235.780 ± 434.839 


0.0109 


-507.897 ± 242.718 


0.0210 


-147.202 ± 76.589 


0.0404 


-27.551 ± 14.666 


0.0777 


-1.677 ± 0.894 


0.1494 


-0.126 ± 0.035 


0.2872 


-0.004 ± 0.001 



(a) Bispectrum of the quantity Sp = Ff < F > —1, where 
F is the observed continuum-fitted flux and < F > is the 
mean flux for each QSO, as computed from eq. (1). 

• on scales k < 0.03 s/km an analytical relation be- 
tween the flux bispectrum and flux power spectrum 
based on second-order perturbation theory in the frame- 
work of the FGPA reproduces the slope but not the am- 
plitude of the bispectrum. 

In summary, the observed bispectrum obtained 
from the LUQAS sample agrees similarly well with that 
of absorption spectra obtained from numerical simula- 
tions and randomized observed absorption spectra. This 
suggests that significantly larger samples of observed 
spectra as for example expected from SDSS and a tight 
control on systematic errors are necessary to utilize the 
bispectrum to constrain the 3D distribution of the ab- 
sorbers and/or cosmological parameters. 
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